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Raman spectroscopy is a useful experimental tool to investigate local deformation and structural 
changes in SiCb-based glasses. Using a semi-classical modelling of Raman spectra in large samples 
of silica glasses, we show in this paper that shear plastic flow affects the Raman measurement in 
the upper part of the spectrum. We relate these changes to structural modifications, as well as a 
detailed analysis of the vibration modes computed in the same frequency range. These results opens 
the door to in situ monitoring of plastic damage in silica-based structures. 
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I. INTRODUCTION 

At present, amorphous materials continue to play a 
role of key materials in many technological areas, ranging 
from Si-based microelectronics to 5D Opt ical Data Stor¬ 
age by Ultrafast Laser Nanostructuring 1 2 . Among them, 
silicate glasses are one of the most important, for exam¬ 
ple for their use as glass fibers, as dielectric insulator, 
or in the building industry. However these materials are 
brittle at large scale, with possibly catastrophic propaga¬ 
tion of cracks, when they are submitted to large stresses. 
It is now well known that crack propagation in silicate 
glasses is initiated around plastic deformations occurring 
at the sub-micrometer scale 3 . The identification of plas¬ 
tic process zones would thus be of great interest in the 
effort for damage reduction in such materials. However, 
in situ measurements of plastic deformation are difficult 
at the micrometer scale. Thus, an heightened interest 
maintains nowadays to understand the structural atomic 
(re)arrangements in these materials. Diffraction probes 
like x-ray and neutron diffraction applied to amorphous 
materials can only partially characterize their atomic 
structure, since these materials lack translational symme¬ 
try. Therefore the vibrational spectroscopy techniques, in 
particular, Raman, Hyper-Raman, and Brillouin scatter¬ 
ings are applied as a complementary class of experimen¬ 
tal techniques. Recent studies have revealed that Raman 
spectroscopy is not only sensitive to the local structure, 
but is a good local densification sensor in silicate^. In 
polymeric systems submitted to large shear deformation 
with shear banding, it was shown that Raman measure¬ 
ments are correlated to the direction of the shear band 5 . 
In disordered materials, vibrational response is very sen¬ 
sitive to the proximity to plastic instability. However, 
shear deformation in silicate glasses takes place at a very 
small scale, and it is not clear wether the related struc¬ 
tural changes would have a measurable Raman signature. 
Moreover, even with a visible signature of the deforma¬ 


tion, it is rather difficult to extract directly structural in¬ 
formations from the experimental data, as this requires 
an accurate theoretical modelling of the material under 
investigation. 

Our study is aimed to investigate the microscopic ba¬ 
sis of a signature of plastic shear in Si0 2 glasses with 
the help of a modelling of Raman spectra coupled with 
classical molecular dynamics simulations. For crystalline 
solids, the first-principles methods have proved very suc¬ 
cessful in reproducing Raman cross-sections 7 . However, 
the application of these approaches to calculate the Ra¬ 
man spectra of disordered materials has long been pre¬ 
vented by the computational effort required for treating 
systems of large system sizes. Furthermore, it was shown 
by one of us^ that long range correlations take place in 
the inhomogeneous vibrational response of glasses such 
that this size should be bigger than 50 A, that corre¬ 
sponds to ^ 8000 atoms, to avoid the size-effecti 9 . The 
paper is organized as follows: section 2 describes the 
preparation of the glass samples; section 3 gives details 
about semi-classical Raman simulations, and the calcu¬ 
lation of the vibrational properties; in section 4 we show 
the results obtained on the calculation of Raman spec¬ 
tra, including a comparison with experimental results, 
and the sensitivity to plastic shear, with a description of 
the related vibration modes and structural changes. The 
last part serves to an additional discussion of the results 
in the light of previous modelling. 


II. GLASS MODEL 

The amorphous silica model used in our numerical 
study was prepared following our previous work 10 . An 
amorphous system of N = 8232 atoms is obtained within 
classical molecular dynamics simulations using a melt- 
and-quench procedure. The model is contained in a box 
with linear size 50.16 Asufficiently large to get acous- 
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tic plane waved^. The equilibration of the liquid, quench 
and relaxation of the glass are performed classically using 
the van Beest, Kramer and van Santen (BKS) potential 11 
modified by A. Carre et alH It can be described as a 
two-body potential 

+ $«p ck (r), (1) 

where a and /? are the types of atoms (O or Si) and r is 
the distance between them; 

$a/T*( r ) = Qaqise 2 V w (r)G w (r), (2) 
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0-0 

1388.773 
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175.0 
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Si-0 

18003.7572 
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Si-Si 
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r inf (A) 


0-0 

-14.97811134 

39.0602602165 

1.75 


Si-0 

-3.24749265 

-15.86902056 

1.27 


Si-Si 

0.0 

0.0 

0.0 



TABLE I: Parameters of the empirical potential used to model 
the silica glass. 
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where j sh = yw = 0.5, r c , w = 10.17 A, r c , sh = 5.5 A. 

We also add a strong and regular repulsive part at short 
range (r < r* n /) to avoid the collapse of atoms at high 
pressure, or high temperatures. The added repulsive part 
has the following form: 


O r ) = (^f ) 12 + E °pr + F a p. (7) 

D aP , E a p , and F a p have been adjusted in order to have 
the continuity of the potential and its first, and second 
derivatives. The parameters of this potential are tabu¬ 
lated in Table HI 

The total energy can then be written 

Etot = E *l*? W (ra) + ^j) ( 8 ) 

i<j 

and will be used to compute the Dynamical Matrix, as 
will be explained later. 

The sample preparation proceeds as follows. The 
glasses are obtained from a crystal state which is heated 
and then quenched. We start from a /3-Cristobalite crys¬ 
tal sample, which is the last crystalline state before the 
melting point at normal pressure. Then we heat the 


sample up to 5200 K during 1 ns, let it evolve at con¬ 
stant temperature during 10 000 time steps (time step 
St = 10“ 15 s), and quench it at 0 K in 1 ns. The quench¬ 
ing rate is hence 5.2xlO 12 K s _1 . Finally we relax the 
simulation box to avoid residual stress, and then obtain 
a density of 2.18 g cm -3 with a pressure about 0 GPa 
with a precision of 10 MPa. The structural properties 
of such system have already been studied extensively 10 . 
We show in Fig. [l] the distribution of Sz-rings in the sys¬ 
tem. The Sz-rings describe the medium-range ordeP 3 of 
vitreous SiC> 2 . They are defined as the smallest closed 
loop obtained starting on a Si atom and following Si-0 
bonds 14 . A 3-membered (or 3-fold) ring contains three Si 
atoms, a n-membered (n-fold) ring contains n Si atoms. 
The distribution of rings in silica glasses has often been 
related to the amplitude of D\ and D2 bands in Raman 
spectrsPKEl. In Fig. [l] we see that for the sample stud¬ 
ied at rest in this paper, the distribution is maximum 
for 6-membered rings, and the proportion of 5- and 7- 
membered rings is approximately the same. This distri¬ 
bution is sensitive to the pressure, the density, and the 
quenching rate used during the sample’s preparation. We 
see in this paper that it depends as well on the plastic 
shear (Fig. [2]), with an increasing number of small rings 
upon shear. It should be noted, that previously^, the 
glass-samples prepared for a number of particles smaller 
than 100, had significant variations in rings distribution 
conditionally on the time spent in the liquid state. How¬ 
ever, it was mentioned that this distribution depends on 
the number of particles 19 . We found out that the vari¬ 
ation becomes weak (less sensitive to the time spent in 
the liquid state during the preparation of the glass (see 
Fig. 0 with an increasing number of particles. The ob¬ 
tained distribution is in agreement with the last results 
of D. Waroquier^Sl, and a similar distribution was also 
obtained in the recent experimental results of monolayer 
SiOj™. For unknown reasons, the number of small 
3- and 4-membered rings is always larger in our systems 
than in experimental measurements. We will comment 
later it’s evolution upon plastic shear deformation. The 
distribution of Si-O-Si angles is shown in Fig. [3] The 
total distribution is in agreement with previous resultd^ 
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with a maximum at 150° db 20° and agrees reasonably 
well with experimental results^HHI The distribution of 
Si-O-Si angles in 3-fold rings is more peaked and cen¬ 
tered on lower values (135° ±5°), while the distribution 
of Si-O-Si angles in 4-fold rings is very close to the total 
one. Note that the values obtained in our samples for 
the distributions restricted to 3-fold and 4-fold rings are 
larger that those obtained in the preceding work of A. 
Pasquarello et alJ 16 and closer to the experimental one. 
Upon shear, the distribution shows a secondary distribu¬ 
tion centered on 100°. This auxiliary angle corresponds 
to adjacent tetrahedra, and was already seen upon com¬ 
pression 10 but in larger proportions. Upon shear, the 
effect is only slightly visible, and affects mainly Si-O-Si 
angles belonging to rings whose size is larger than 4-fold. 
In the same time, it was shown (see Fig. [5| that only a 
few percent of particles vary its coordination number: 3% 
additional 5-coordinated Si-atoms (as discussed id22) ap¬ 
pear upon shear in the sample. Local structural changes 
are thus only tiny visible upon shear at constant volume. 
We will discuss later the more pronounced evolution of 
mesoscopic structural correlations (evolution in the num¬ 
ber of rings) upon shear. 

We have chosen to perform the mechanical deforma¬ 
tion of our systems in the a-thermal regime, i.e. we ei¬ 
ther perform energy minimization (T = 0 K), or we let 
our glass evolve at a very low temperature T = 10 -5 K, 
temperature at which thermal effects are totally negligi¬ 
ble and not high enough to activate plastic events^. For 
the shear at constant volume, we changed the shape of 
the box within a monoclinic symmetry and a small defor¬ 
mation step, and then let the atoms relax by minimizing 
the total energy in the box at constant volume. We used 
the Polak-Ribiere version of the conjugate gradient (CG) 
minimization algorithm proposed by LAMMPS. At each 
iteration, the force gradient is combined with the pre¬ 
vious iteration information to compute a new search di¬ 
rection perpendicular (conjugate) to the previous search 
direction. We carry on this operation to get the desired 
total shear deformation. The shear step corresponds to 
a deformation of Ae xy = 10 -4 , that is small enough 
to decorrelate the plastic events that occur during the 
shea]! 27 l It corresponds to a very fast relaxation time 
between each strain step (quasi-static deformation) that 
allows to follow a mechanical equilibrium path (athermal 
response). This kind of deformation gives very detailed 
results on the restructuring of the phase space upon me¬ 
chanical deformatiorP 6 !, and allows to construct the land¬ 
scape over which thermal activation acts. Thermal acti¬ 
vation, that can play a role even below the glass transi¬ 
tion temperature^, will thus not be discussed here. 

The resulting stress-strain curve is shown in Fig. [4] It 
shows a linear response followed by stress-softening and 
hatched stress relaxation processes. At large deformation 
(> 30%) the system reaches a visco-plastic plateau char¬ 
acterized by a constant average visco-plastic flow stress. 
The hatches, that are the fluctuations around the aver¬ 
age flow, correspond to local plastic rearrangements giv- 
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FIG. 1: Ring statistics for uncompressed a-SiCL sample 

(8232 particles). 
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FIG. 2: Rings statistics in a-Si02 upon shear. Comparison 
between the initial distribution of rings (red) and the distri¬ 
bution of rings after the sample entered in the plastic plateau 
(black). 


ing rise to energy dissipation 10 . We compared the results 
obtained on 10 initial configurations at rest, and 10 con¬ 
figurations chosen randomly in the plastic plateau. 

III. SIMULATION DETAILS 
A. Raman model 

Here we briefly outline the formulation that we used 
for the calculation of Raman activities. We only focus 
on first-order processes, which involve a single phonon 
excitation. For a model system, the vibrational frequen¬ 
cies uo n and their associated eigenmodes are found by 
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FIG. 3: Distribution of Si-O-Si angles (a) in 3-fold rings, (b) in 4-fold rings, and (c) total distribution. 



FIG. 4: Shear stress as a function of engineering shear strain 
in a-Si 02 sample upon quasi-static and athermal shear defor¬ 
mation. 


solving the set of linear equations: 

(9) 

Bj 

where D is the dynamical matrix, and where the capital 
Latin indices run over the atoms, and the uppercase Latin 
indices are the three Cartesian directions, respectively. 
The dynamical matrix D is defined by 

t r) _ 1 _ d 2 E tot 

Al ’ Bj V 7 Mjm dR A idR B j ’ 1 ’ 

where E tot is the total energy of the system, and Ra 
and Ma are the position and the atomic mass of atom A. 
Diagonalization of the dynamical matrix D gives the vi¬ 
brational frequencies cj 2 and the normalized eigenmodes 
£ai- The Raman susceptibility tensors 28 7are then 
given for each mode n by 



strain, % 


FIG. 5: Evolution of the coordination number, and rings 
distribution in a-Si 02 model sample, upon athermal plastic 
shear deformation. 
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dxii q fc 

BRah 


( 11 ) 


where x 1S the electric polarizability tensor and V is the 
volume of the cell. 

In experimental set-ups, it is customary to record the 
Raman spectra in the HH and the HV configurations, 
in which the polarization of the outgoing photons is re¬ 
spectively parallel and orthogonal to the ingoing photon 
polarizatioEpS. Using the isotropy of disordered solids, 
we express the contributions of the n-th mode to the HH 
and HV Raman spectra aP^ 


T n — n 2 + 
l HH — a n + 

( 12 ) 

jn _ ^ 7,2 

1 HV - 

( 13 ) 


where a n and b n are obtained from 
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d n 


3 (^ 11 + 7 * 22 +^ 33 ). 


( 14 ) 


^ = ^[(^ 11 -^ 22 ) 


+ (^ 22 -^ 33 ) 


+ 3 


+ (^11 - ^ 33 ) + 

(^ 2 ) 2 + ( 7^ 3 ) 2 + ( 7^ 3 ) 2 


( 15 ) 


Focusing on the Stokes process, in which a vibrational 
excitation is created by an incoming photon, we express 
the total power cross-section as (in esu units) 


where R is a unit vector along R = Ra — Rb, and a[ 
are the derivatives of the bond polarizabilities with re¬ 
spect to the bond length (a[ = (doq :P /dR)\R = R 0 , where 

Ro is a typical distance). Therefore, when one type of 
bond occurs, the bond polarizability model is completely 
defined by three parameters: 2a p + <aj, a[ — a p , and 
(cq — a p )/R. We used the parameters of the bond polar¬ 
izability model already derived by one of us in R e f. 3 2 iMl, 
whose values are summarized in Table ED 


TABLE II: Bond-polarizability-model parameters. 


T P ( \ a ( Ul ~ ^n) 4 ^ 

I HH,Hvi u ) = 47r - q4 - X 

n 

x IhH.HVT, + 1]^(^ — ^n)> (16) 

ZUJn 

where V is the volume of the scattering sample, ujl the 
frequency of the incoming photon, c the speed of light, 
and n{(jj) the boson factor: 

n(uj) = -——^——-. (17) 

exp(huj/ksT) — 1 

The bond polarizability model 30 31 has successfully 
been applied for the calculation of Raman intensities in 
a large variety of system^. In this approach, the polar¬ 
izability is modeled in terms of bond contributions: 


Model 2a p + cq a[ — a' p (cq — a p )/R 


BP 

0.771 

0.196 

0.056 

((47r) _1 -Bohr _1 ) 





B. Vibrational density of state 

The frequencies uj n and the corresponding normalized 
eigenmodes ^ are obtained by diagonalizing the dy¬ 
namical matrix. The Intel MKL tool was used for the 
diagonalizatiorH The associated atomic displacements 
are given by: 


XiM) = hY, a v( A ’B) (18) 

B 

with 


(7l, B) — 2 (2 Qtp + OLi)8ij-\- 

+ («/-«p)(^-^i), (19) 

where R = Ra — Rb is a vector which defines the di¬ 
rection and the distance of a pair of nearest neighbor 
atoms at sites Ra and Rb The parameters cq and a p 
correspond to the longitudinal and perpendicular bond 
polarizability, respectively. 

The bond polarizability model further assumes that 
the bond polarizabilities cq and a p only depend on the 
length of the bond. Thus, the derivative of the local bond 
polarizability with respect to the relative displacement of 
the atoms A and B yields: 


daij{A,B) 

dRAk 
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(2a' p + 5ijRk~\~ 
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( 21 ) 


The index n labeling the vibrational modes runs from 
1 to 3V, where N is the total number of atoms in the 
model. 

The u-DOS can be decomposed according to the 
weights of the two species: Z{uS) = Z a (uS). The par¬ 
tial density of states Z a (u )) is defined by: 


N a 3N 

= ^(u-COn), (22) 

A n 

where A runs over the set of atoms belonging to a given 
species a. The u-DOS can be decomposed as well ac¬ 
cording to the weight of different kinds of vibrations 
(e.g. bending, rocking, stretching). The vibrations of 
the bridging oxygens (BO) may be classified as stretch¬ 
ing, bending, or rocking according to whether the BO 
atoms move parallel to the line joining the the two Si 
neighbors, parallel to bisector of the Si-O-Si angle or per¬ 
pendicular to the Si-O-Si plane, respectively. Thus, if we 
apply following unit vectors: r 0 the unit vector perpen¬ 
dicular to the local Si-O-Si plane, fi + i/, the unit vector 
parallel to the bisector of the Si-O-Si angle, and n _v 
the unit vector perpendicular to r 0 and fi + i/, then the 
corresponding projection of the vibration modes in the 



















rocking, stretching and bending and contributions con¬ 
tribute to the partial density of states as 
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Zo(u>) — ^rockiyj') H - ^stri^) H - ^ beni (23) 
with a retricted to Oxygen atoms, and 


1 No 3N 

^rocfeH = — Xll^A-r©| 2 (5(w-w n ), (24) 

A=Ox n 
1 No 3N 

Zstr{u) = 5(id — LU n )i (25) 

A=0 x n 
1 N 0 3N 

Zbe„(u) = — ^2 H^A-ri+i'| 2 (5(w-w n ). (26) 

A=Ox n 


IV. RESULTS 

A. Validity of the semi-classical calculation 

We obtain the Raman HH and HV spectra by com¬ 
bining the classical eigenmodes obtained with the diago- 
nalization of the (large scale) dynamical matrix and the 
BP model to describe the vibration dependence of the 
local polarizability. We computed the Raman spectra 
for the uncompressed samples (see Fig. [6] and [7]) and in 
the pressure imposed case (see Fig. [8]) to validate our 
model. The results are satisfying. Note that the spectra 
obtained in these figures have been obtained on a single 
sample without ensemble average. The computed spec¬ 
trum was simply smoothened with gaussian functions of 
width 20 cm -1 close to the experimental resolution. We 
see in Figures |6j[7| and [8] that the main characteristics of 
the experimental spectra are recovered within our semi- 
classical approximation. 

In the parallel polarization case (HH spectrum), we 
can recognize well the main band (^ 400 - 550 cm -1 ), 
and we can see, that the high-frequency experimen¬ 
tal bands, in particular those located at ^ 800 cm -1 , 
~ 1050 cm -1 , and ~ 1200 cm -1 are well reproduced. 
The cross-polarized HV Raman spectrum gives similar 
results. These results are in good agreement with the ex¬ 
perimental measurement^, despite only the small shift 
of the 800 cm -1 band to the left in the HV spectrum, 
with respect to the experimental data. In the pressure 
imposed case, the width and the position of the main 
band (400 - 550 cm -1 ) has the same trend as in the ex¬ 
perimental spectra (see Fig.[8|: it shows a shift to higher 
frequencies and a progressive narrowing, that were the 
main results oP^l This successful comparison with ex¬ 
perimental data encourages us to perform further anal¬ 
ysis of our semi-classically computed Raman spectra for 
simple shear deformation of glasses, despite the lack of 
experimental data in this case. 



FIG. 6: Raman HH spectrum in uncompressed a-SiCb- Red 
(solid) curve is theory [current work], black (dash) curve is 
experimental data from 35 . 



FIG. 7: Raman HV spectrum in uncompressed a-Si02. Red 
(solid) curve is theory [current work], black (dash) curve is 
experimental data from 35 . 


B. Shear sensitivity of Raman spectrum 

We mentioned above that the quasi-static deformation 
of our system is imposed through very fast energy relax¬ 
ation after small imposed shear strain steps. The cor¬ 
responding stress-strain behaviour is depicted in Fig. [4j 
We are interested in the plastic regime, and more pre¬ 
cisely in the plastic flow depicted in the plastic plateau. 
It seems reasonable to average over the configurations 
from the plateau of the stress-strain curve, since it is dif¬ 
ficult to say what is the exact stage of our considered 
sample upon shear at constant stress in indentation ex¬ 
periment. The corresponding Raman HH and HV spec¬ 
tra, averaged over 10 configurations, are shown in Fig. [9] 
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FIG. 8: Raman HH spectrum in a-Si02 after an isotropic de¬ 
formation has been applied, leading to different internal pres¬ 
sures P=0 GPa, P=3 GPa, and P=5 GPa (8232 particles), 
and comparison with experimental measurements from--. 
The curves have been shifted vertically by an arbitrary value 
for clarity. 



FIG. 9: Raman HH spectrum in a-Si02 upon shear. Black 
(solid) - basic sample, red (dash-dot) - mean spectrum, blue 
(doted area) - different configurations, orange (dot) - fluctu¬ 
ations. 


and [10} respectively. Blue-dots is a coverage from con¬ 
tributions of different spectra which we average over, red 
(dash-dot) curve is the mean spectrum obtained on the 
sheared samples, and the black line is the spectrum of the 
initial unsheared sample. We see that the fluctuations 
(dotted curve) are not to large, that is an indisputable 
advantage of using larger samples in current work (8232 
particles) with respect to most past works, where smaller 
(~ 1000 particles) samples could be analyzed hardly due 
to finite size effect 19 . In our semi-classical simulations of 
Raman spectra performed on sufficiently large systems, 
it appears easy to separate the effect of shear from the 



FIG. 10: Raman HV spectrum in a-Si02 upon shear. Black 
(solid) - basic sample, red (dash-dot) - mean spectrum, blue 
(doted area) - configurations, orange (dot) - fluctuations. 


configurational noise. 


Both HH and HV spectra show clearly the activation 
of high-frequency modes upon shear, whereas the low- 
frequency domain (below 500 cm -1 ) remains without 
valuable changes. Both spectra show an increase of the 
Raman activity in the range 800-1000 cm -1 as does the 
density of vibrational states (Fig. 12), and a shift of the 
1050 cm -1 band to lower frequencies. An analysis of 
the high-frequency mode contributions is not obvioud^, 
therefore we must be careful to treat this range. Our 
simulations show that the plastic shear flow induces irre¬ 
versible structural changes in the model silica glass stud¬ 
ied here, with a transformation of large rings into lower 
order units. The concomitant increase of the D2-band 
could support the usual interpretation in terms of vibra¬ 
tions of 3-fold ring^l 


At larger frequencies, the range 800 - 1000 cm 1 has 
an obvious growth for HV & HH spectra. Particularly, 
the HV spectrum shows clearly this growth, whereas the 
bands ^ 1050 cm -1 and ^ 1200 cm -1 become less ex¬ 
pressed with shear but their contribution to the spectra 
remains. Note that the Raman signal is not sufficiently 
sensitive to discriminate between the configurations of 
silica in the plastic plateau that are close to a plastic 
instability, and those that are deeply stable, unlike irP. 
This is probably due to the fact that Raman signals are 
more sensitive to the irreversible structural changes that 
appeared in the beginning of the plastic plateau, than 
to the sparse mechanical instabilities responsible for the 
small scale stress fluctuations in the plateau. We will 
now relate this behaviour to a detailed study of the shear 
sensitivity in the density of vibrational states, and then 
to a detailed study of the vibration modes at different 
eigenfrequencies. 
























FIG. 11: VDOS spectrum in a-Si02 upon shear. Black 

(solid) - basic sample, red (dash-dot) - mean spectrum, blue 
(doted area) - configurations, orange (dot) - fluctuation, cross 
- experimental results from 37 . 



FIG. 12: Basic vs Shear of VDOS spectra in a-Si02 with the 
different contributions of Si and O atoms to the vibrations. 
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FIG. 13: Basic vs Shear of VDOS spectra in a-Si02 with 
the different contributions of bending, stretching, and rock¬ 
ing types of vibrations. The total corresponds to the oxygen 
contribution shown in Fig. [12] 



C. Density of Vibrational States 


FIG. 14: Basic vs Shear of VDOS spectra in a-Si02 with the 
different contributions of 3-fold and 4-fold rings. 


Following the same logic as before for Raman, the Vi¬ 
brational Density of States (VDOS) can be evaluated 
with the help of Eq. (22). We can see in Fig. 11 the 


changes appearing in the VDOS after the system entered 
in the plastic flow regime. These changes are mainly 
located in the high-frequency domain (> 1000 cm -1 ) 
and in the intermediate frequency range (700 cm -1 - 
1000 cm -1 ), that correlate with the Raman spectra. Due 
to the activation of the VDOS in the vicinity of ^ 950 
upon shear, the populations of Raman HH and Ra¬ 


cm 


-l 


man HV spectra have the same trend. The changes in the 


VDOS are related to a redistribution of the vibrational 
activity: on one hand we observe the decrease of popula¬ 
tion in the ^ 500 and ^ 1000-1200 cm -1 bands, on the 
other hand an activation of the mute area between these 
bands occurs, as well as an enhancement in the range of 
D2-band (^ 600 cm -1 ). 

The analysis of the partial contributions from Si-/0- 
atoms to VDOS (see Fig. 12) is very instructive . In 
agreement with previous results obtained in Refs.®^ 
on different simulation of silica samples, it shows that 
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FIG. 15: Basic vs Shear of VDOS spectra in a-SiC >2 with the 
different contributions of Si and 0 atoms to the vibrations 
supported by 3-fold rings. 


FIG. 17: Basic vs Shear of VDOS spectra in a-Si02 with the 
different contributions of Si and O atoms to the vibrations 
supported by 4-fold rings. 



FIG. 16: Basic vs Shear of VDOS spectra in a-Si02 with 
the different contributions of bending, stretching, and rocking 
types of oxygen atom vibrations supported by 3-fold rings. 



FIG. 18: Basic vs Shear of VDOS spectra in a-Si02 with 
the different contributions of bending, stretching, and rocking 
types of oxygen atom vibrations supported by 4-fold rings. 


the high frequency (> 950 cm -1 ) and the low frequency 
domains (< 650 cm -1 ) for the samples at rest, are dom¬ 
inated by the vibrations of oxygen atoms. A peak in the 
contribution of Si displacements to the VDOS is located 
at ^730 cm -1 . Upon shear, the decay of the high fre¬ 
quency VDOS is univocally related to the decay in the 
oxygen motion, with a poorer definition and a small shift 
of the two peaks located at 1050 cm—1 and 1200 cm -1 . 
Similarly, the increase of the activity in the intermediate 
range 700-1000 cm -1 is due mainly to the increase of the 


oxygen displacements in that frequency range. 

Finally, the decomposition of the modes into bend¬ 
ing/stretching/rocking modes (see Fig. 13 ) follo wing the 
Eqs. ( [24] ) - ( |26| ) shows in agreement wrt n 13|16|39|4Q | that 
rocking and bending modes dominate at low frequencies, 
while the high frequency range is dominated by stretch¬ 
ing modes. Upon shear, we show in this figure, that the 
increase of oxygen activity in the intermediate domain 
is uniquely due to the increasing number of stretching 
modes in that frequency range 
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FIG. 19: Basic vs Shear contribution to the VDOS of in- 
phase breathing modes supported by oxygen atoms in 3-fold 
and 4-fold rings, divided by the corresponding number of 
rings. 


In order to compare with the usual interpretation of 
high-frequency vibrations, we will now study in more de¬ 
tails the contribution of the different rings to the vibra¬ 
tions, and its dependence with the plastic shear deforma¬ 
tion. Fig. p| shows the contribution of 3-fold and 4-fold 
rings to the VDOS. The general increase of 3-fold and 
4-fold contribution upon shear, is due to the increasing 
number of such rings upon plastic shear (Fig. [2|. In the 
high frequency domain, the contribution of 4-fold rings 
is stable, but the increase of the contribution of 3-fold 
ring is opposite to the decay of the total VDOS in this 
frequency range. This proves that the decay of the oxy¬ 
gen displacements responsible for the decay of the total 
activity at high frequency is located outside the small 
rings. 

Let’s first focus on the contribution of 3-fold rings 
to the activity. Fig. [15] shows the partial contribution 
of 3-fold rings to VDOS, normalized by the number of 
rings. The main contribution of 3-fold rings appears at 
^770 cm -1 (thus at a frequency higher than the fre¬ 
quency of the D2 band). This peak is mainly due to 
the contribution of bending modes (Fig. [Tg| ) and it is 
smoothened upon shear, for oxygen as well as for sili¬ 
con atoms. Oxygen atoms show however an increasing 
activity in the intermediate frequency range (700 - 1000 
cm -1 ). This increase of the activity in the intermediate 
frequency range is related to an increase of stretching 
modes, as already seen for the whole system. 

Let’s consider now 4-fold rings. The main contribu¬ 
tion of 4-fold rings is located at ^700 cm -1 and is due 
tot the enhanced vibration of silicon atoms at this fre¬ 
quency (Fig. ED- This peak is smoothened upon shear, 
but there is a concomitant increase of the oxygen activity 
in the intermediate frequency range (700 - 1000 cm -1 ) 
upon stretching modes, that contributes to the general 
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FIG. 20: Basic vs Shear Participation Ratio in a-Si02, as 
a function of the frequency, and partial Participation Ratios 
restricted to O and Si atoms. 
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FIG. 21: (a) Vibration mode at the frequency 97 cm -1 for 

the basic configuration and (b) at the frequency 849 cm -1 . 
The arrows are proportional to the atomic displacements. Red 
arrows are for oxygen atoms, and Black arrows for silicon 
atoms. 


increase of VDOS in the intermediate frequency range. 
Note interestingly, that the main contributions of oxy¬ 
gen are located at 150 cm -1 and at 630 cm -1 , and do 
not correspond to D1 band. 

In order finally to test the interpretation given by 
GaleeneJ 15 -^ for D1 and D2 band in our samples, we 
have isolated the contribution of in-phase bending vibra¬ 
tions of oxygen atoms in 3-fold and 4-fold rings, that 
correspond to the breathing vibration modes of the rings. 
The contribution of these modes to the VDOS is shown 
in Fig. [19] It is remarkable that we see two well defined 
peaks. It should be noted that the in-phase breathing 
modes of Oxygen atoms in 4- and 3-fold rings have been 
attributed respectively to the defect lines D1 and D2 in 
the HH Raman spectrum at 495 and 605 cm -1 . It is not 
a surprise that we find the two lines at higher frequen¬ 
cies as similar overestimations 41 are found when the BKS 
potential is used for evaluating vibrational modes of sil¬ 
ica polymorphs for which experimental data are known. 
However, we show here that we can predict the behaviour 
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FIG. 22: Contribution of Si vs O atoms to the displacement 
field in the eigenvectors, as a function of the mode’s eigenfre- 
quency. Comparison upon steady shear. 



FIG. 23: Contribution of 3- & 4-fold rings to the displace¬ 
ment field in the eigenvectors as a function of the mode’s 
eigenfrequency, and effect of shear. 


of D1 and D2 lines in the Raman HH spectrum upon an 
external perturbation such as shear. Indeed, D1 and D2 
lines are expected to increase upon shear as can be seen 
in Fig. [9] keeping in mind the we assign the theoretical 
D1 line to the feature at ^630 cm -1 and the D2 line to 
that at ^770 cm -1 . 


D. Analysis of the Vibration Modes 


The vibration modes are the resonant eigenmodes ob¬ 
tained by the exact diagonalization of the Dynamical Ma¬ 
trix, as explained in Part |nTAl 

A parameter that is frequently used to characterize the 
vibration modes in amorphous systems is the Participa¬ 
tion Ratio (PR). It is defined as 


PR = 


1 (Ef = 1 |u ,| 2 ) 2 

N £f = 1 kl 4 


(27) 


It quantifies the percentage of particles moving with 
the same amplitude (PR = 1/N if one single particle 
is moving, PR = 1 if the N particles are moving to¬ 
gether), and is often used to characterize the degree of 
localization of the modes. The Participation Ratios of 
the vibration modes obtained from the eigenvectors of the 
dynamical matrix are depicted in Fig. [20] It shows dif¬ 
ferent frequency regimes: the first 3 translational modes 
in each of the 3 directions have a Participation Ratio, 
PR = 1. At very low frequencies (< 50 cm -1 ), the Par¬ 
ticipation Ratio of the plane waves decays with the fre¬ 
quency (smaller wavelength). In the same low frequency 
regime (< 50 cm -1 ), we can identify soft modes that are 
characterized by a very low Participation Ratio. These 
modes correspond to localized vibrations with unusual 


low frequencies, that are precursors of local plastic in¬ 
stabilities. Surprisingly enough, in silica glasses, these 
modes are present whatever the distance to a plastic in¬ 
stability, unlike in other glassy systems where they ap¬ 
pear only at the onset of plasticity. We will not study 
them in more details here, because they are not in the fre¬ 
quency range where the main changes to the Raman spec¬ 
tra appear upon plastic deformation. At higher frequen¬ 
cies, close to 100 cm -1 , the Participation Ratio reaches 
a local maximum, followed by a plateau. The frequency 
where this local maximum appears (with PR « 0.45) is 
close to the Boson’s peak frequency obtained by divid¬ 
ing the VDOS by the Debye predictiorl 10 . At this fre¬ 
quency, the modes display large scale rotational motions 
structures (see Fig. [2l]-a). In the plateau, the rotational 
motion is more noisy. At higher frequencies, the Partici¬ 
pation Ratio shows a small decay and a local maximum 
(uj ^850 cm -1 ) followed by a strong decay at a frequency 
close to 900 cm -1 . The global Participation Ratio dis¬ 
cussed here, can be compared to the partial Participation 
Ratios, where the sum in Eq. (27) runs only over one type 
of atoms (Si or O). It is seen in Fig. 20, that the Participa¬ 
tion Ratio restricted to oxygen atoms is very close to the 
global one. The Participation Ratio of silicon atoms does 
not depend to much on the frequency and is quite high, 
indicating homogeneous motion of the silicon skeleton. 
The variation of the Participation Ratio near 850 cm -1 
upon shear is directly related to a change in the spa¬ 
tial distribution of oxygen displacements, that start to 
localize smoothly at 700 cm -1 while the silicon matrix 
still vibrates collectively (PRsi ^0.55 but PRo ^0.25 
at 800 cm -1 ). Oxygen motion appears thus very sensi¬ 
tive to the frequency in the intermediate frequency range 
700-900 cm -1 . After a strong decay of the Participa¬ 
tion Ratio corresponding to localization processed^, the 
PR increases again between 1000 and 1200 cm -1 . As 
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FIG. 24: Separate contribution of Si & O atoms in the 3-fold 
rings, to the displacement field in the eigenvectors as a func¬ 
tion of the mode’s eigenfrequency, and effect of shear. The 
curves corresponding to the sheared state have been normal¬ 
ized by the ratio between the number of 3-fold rings in the 
sheared state, and the number of 3-fold rings in the initial 
state. 


FIG. 25: Separate contribution of Si & O atoms in the 4-fold 
rings, to the displacement field in the eigenvectors as a func¬ 
tion of the mode’s eigenfrequency, and effect of shear. The 
curves corresponding to the sheared state have been normal¬ 
ized by the ratio between the number of 4-fold rings in the 
sheared state, and the number of 4-fold rings in the initial 
state. 


seen before, these high frequencies correspond to local 
stretching modes 43 . The main effect of plastic shear is 
to smoothen strongly the bump in the participation ratio 
close to 850 cm -1 . This effect is related to a partial local¬ 
ization of the modes supported by oxygen atoms in this 
frequency range upon plastic shear, and will be discussed 
later in more details. 


In order to describe more precisely the vibration modes 
in the different characteristic parts of the spectrum, we 
have decomposed the vibration modes into different con¬ 
tributions: the displacements supported by the different 
types of atoms (O and Si, Fig. [22]) , the displacements sup¬ 
ported by different kinds of rings (3-fold or 4-fold rings 
in Fig. 23), and the contribution of the vibrations of the 
different atoms O and Si in the different rings (Figs. [24] 
to 27). In Fig. 22, we show the amplitude of the dis¬ 


placements supported by the different types of atoms in 
each vibration modes as a function of the frequency. We 
see clearly in this figure, that the contributions of the 
different atoms to the vibrations change mainly in the 
700-900 cm -1 range, upon plastic shear. Moreover, it is 
seen in this figure, that there is a transfer from the dis¬ 
placements initially supported by Si atoms that freeze, 
to displacements of O atoms that increase. 


We will now look at the contribution of the vibrations 
in the different Si-rings. In the usual interpretations of 
Raman spectr#* 5 F 6 l 44 l 45 ( D1 band is related to vibrations 
of 4-fold rings, and D2 band to the breathing mode of 3- 
fold rings supported by Oxygen atoms, assuming that 
Si atoms are fixed. However, it should be noted, that 
to assign definitely the signatures of 3-fold and 4-fold 
rings is still questionable 19 . This was initially related to 


calculation of eigenfrequencies of isolated rings®®!, and 
then obtained_with ab initio calculations performed on 
small system^ 16 ! In Fig. 23, we see that the vibrations 
of 4-fold rings support approximately half of the total 
displacement, and that 3-fold and 4-fold rings contribute 
to the displacement field in all the modes and at all fre¬ 
quencies. This can be simply understood by the fact 
that, in amorphous materials, the vibration modes are 
in general collective modes, due to the mechanical homo¬ 
geneity that takes place already at the nanometer scale 46 . 
Even modes with low mean-free path are referred to as 
” quasi-local vibration modes” consisting in the superpo¬ 
sition of rotational waves. Despite the case of isolated 
’’soft modes” occurring at very low frequencies, it is in 
general very difficult to relate a vibration mode to the 
vibration of an isolated entity. Note that the large fluc¬ 
tuations seen in Fig. 23 at frequencies close to 950 cm -1 
and 1250 cm -1 are due to the very poor statistics in that 
frequency range were the number of eigenmodes is very 
low, and to the localization of the vibrations that can - 
or not - be close to a n-fold ring. In Fig. [23] the global 
increase of the curves upon shear is due to the increase 
in the number of the small-size rings. Indeed, the ratio 
between the contribution of the displacements supported 
by atoms in 3-fold and in 4-fold rings is given mainly 
by the ratio between the different number of rings (and 
thus the different number of atoms taken into account 
when computing the partial displacements amplitude). 
A rapid comparison between Fig. [24] and Fig. [25] shows 
that the total displacements supported by silicon atoms 
for example (normalized by the number of rings) is very 
similar for both kinds of rings: it shows a maximum vari- 
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FIG. 26: Separate contribution of Si & O atoms in the 5-fold 
rings, to the displacement field in the eigenvectors as a func¬ 
tion of the mode’s eigenfrequency, and effect of shear. The 
curves corresponding to the sheared state have been normal¬ 
ized by the ratio between the number of 5-fold rings in the 
sheared state, and the number of 5-fold rings in the initial 
state. 


FIG. 27: Separate contribution of Si & O atoms in the 6-fold 
rings, to the displacement field in the eigenvectors as a func¬ 
tion of the mode’s eigenfrequency, and effect of shear. The 
curves corresponding to the sheared state have been normal¬ 
ized by the ratio between the number of 6-fold rings in the 
sheared state, and the number of 6-fold rings in the initial 
state. 


ation upon shear in the 700-1000 cm -1 frequency range. 
When looking to Fig. 23, it is a not surprising to see that 
3-fold rings have a maximum contribution at 800 cm -1 
corresponding to the calculated D2 frequency, while two 
peaks appear for the 4-fold rings: a peak at 175 cm -1 and 
a peak at 700 cm -1 . This last peak should correspond to 
the theoretical D1 band. When looking more precisely to 
the contribution of the different atoms to the vibrations, 
we see a finer structure for oxygen atoms. In 4-fold rings 
(Fig. [25]) , the contribution to the overall displacement of 
the oxygen atoms shows a first peak at 175 cm -1 (ex¬ 
plaining the first peak of 4-fold rings displacements) and 
a secondary peak at 630 cm -1 , close to the theoretical D1 
band, in agreement with the previous analysis (Fig. [I7[ ) 
of the partial density of states. This finer description of 
the vibrations has a signature in the Raman spectrum 
(Fig. [ 9 ]) at the corresponding frequencies. 


In general, the curves showing the contribution of the 
different rings to the displacement field of the eigenmodes 
flatten upon plastic shear. This does not mean that noth¬ 
ing happens inside the sample. The Figures [24] to |27| show 
the separate contributions of Si and O atoms to the dis¬ 
placement field supported by rings of different sizes, and 
the corresponding effect of shear. The curves for the 
sheared configurations have been normalized by the ra¬ 
tio of the number of rings in the sheared vs the number of 
rings in basic configurations, in order to enhance the dif¬ 
ferences only due to shear in the shape of these curves, 
and not due to the quantitative change in the number 
of units. The common feature in all these curves is the 
variation occurring in the 700 cm -1 -1000 cm -1 range. In 
this frequency range, the displacements supported by Si 


atoms decrease, while the displacements supported by O 
atoms increase. This means that upon irreversible plastic 
shear strain, there is an enhancement of the contribution 
of oxygen motion to the vibrations, while the skeleton of 
the rings, supported by Si atoms freezes. This is in agree¬ 
ment with the general trend seen in Fig. [l2j that shows 
that, it the same range, the contribution of oxygen atoms 
to the vibrations increases strongly upon shear, and to 
the localization of vibrations on oxygen atoms seen in 
Fig. [20] The effect of the localization at ^850 cm—1 is 
particularly visible (see Fig. 24) in the oxygen motion in 
3-fold rings. 


V. DISCUSSION 

It was shown in this paper that Raman spectra ob¬ 
tained within a semi-classical simulation of large empir¬ 
ical S 1 O 2 samples (.N > 8 000 atoms) using the Bond 
Polarizability model combined with classical calculation 
of the vibration modes, is able to reproduce the main fea¬ 
tures of experimental Raman spectra, and shows clearly 
a sensitivity of the Raman spectra to plastic irreversible 
shear. The variations of the Raman spectra upon plastic 
shear are located mainly in the intermediate frequency 
range ( 700-1000 cm -1 ) and in the high frequency range 
(> 1000 cm -1 ). In the intermediate frequency range, 
this variation corresponds to enhanced stretching vibra¬ 
tions of oxygen atoms. At smaller frequencies, the global 
increase of the Raman spectrum can be related to the 
increase in the number of small rings upon plastic shear. 
At large frequencies (> 1000 cm -1 ), the contribution of 
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stretching modes supported by oxygen atoms belonging 
mainly to large order rings decays. 

An additional feature concerns the localization of the 
breathing vibrations of oxygen atoms on the 3-fold rings 
at ^800 cm -1 . This unexpectedly large value for the fre¬ 
quency of the 3-fold breathing modes can be attributed 
to an excess of rigidity of the Si-0 bond in the BKS inter¬ 
actions used 41 ! It is however surprisingly enough that the 
D2 and the 800 cm -1 bands (observed experimentally) 
are well reproduced in our approximate calculation of the 
Raman spectrum, and even more surprisingly deserve an¬ 
other microscopic description than the usual one 15 . This 
remark encourages to perform further studies to check 
the validity of the bond polarizability model with a sys¬ 
tematic comparison of the different empirical descriptions 
of interatomic interactions, that is far beyond the scope 
of the present paper. 

The computed Raman spectra are sensitive to plas¬ 
tic shear as soon as the sample entered in the plastic 
flow regime with valuable structural changes. It is not 
sensitive to the first stages of the plastic deformation as 
described by plastic instabilities and soft mode^ whose 
signature affects probably more specifically low frequency 
vibration modes measured by Brillouin scattering. The 
sensitivity to permanent plastic shear may be due to mea¬ 
surable structural changes, especially the increase in the 
number of low-order rings at a mesoscopic scale. This be¬ 
haviour is probably composition dependent and it would 
be interesting to track his signature for other compounds 
like a-Si, silicate glasses, or even Ta^Os glasses known 
to have a very low mechanical noise 47 . Raman spectra 
are however not sensitive to the amount of plastic defor¬ 


mation, unlike in densification experiments. The differ¬ 
ence between shear and densification, is that in the shear 
case, the nature of defects reaches a stationary states, 
and second the average stress applied does not vary any¬ 
more upon plastic flow. Raman spectra can thus be seen 
as stress-sensitive sensors. 

An interesting contribution of semi-classical models, is 
that they allow to connect Raman spectra to a complete 
and detailed description of the large scale vibrational 
modes. This was done explicitly in our paper. Further 
work should concern the application of more accurate 
schemes for assessing vibrational modes. The accuracy of 
them should be first establish through comparison with 
experimental data for systems whose structure is well 
known 41 . This will open the way to very accurate simu¬ 
lation of Raman spectra including position and intensity 
of D1 and D2 defect lines. Moreover, perspectives to this 
work includes in situ measurement of compressive/shear 
plasticity at the micrometer scale in silica glasses, that 
are of crucial interest for crack prediction. 
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